## Plots paper ####

library(tidyverse)
library(ggtext)
library(ggpubr)

ven <- readRDS("ven_elec_2006_2024_final.rds")

vote_share <- ven %>%
  group_by(year) %>%
  summarise(
    p_chavismo = weighted.mean(of_p, validos, na.rm = TRUE),
    p_oposicion = weighted.mean(op_p, validos, na.rm = TRUE),
    p_others = weighted.mean(otro_p, validos, na.rm = TRUE)
  ) %>%
  pivot_longer(cols = c(p_chavismo, p_oposicion,p_others), 
               names_to = "Political_Group", 
               values_to = "Support") %>% 
  # Plot both trends
  ggplot(., aes(x = year, y = Support, color = Political_Group, group = Political_Group)) +
  geom_line(linewidth = 1) +
  geom_point() +
  geom_hline(yintercept = 50, linetype = "dashed", color = "black") +
  scale_color_manual(values = c("p_chavismo" = "red", "p_oposicion" = "blue", "p_others" = "purple"),
                     labels = c("Chavismo", "Opposition", "Others")) +
  labs(title = "Figure 3: Share of the vote by political bloc",
       subtitle = "Venezuela Presidential elections 2006-2024",
       x = "",
       y = "Vote (%)",
       color = "Political bloc") +
  theme_bw() + theme(axis.title.y = element_text(vjust = 2, size = 12),
                          axis.text.x = element_text(size = 10),
                          legend.text = element_text(size = 10),
                          title = element_text(size = 14))
vote_share

ggsave(vote_share, file = "year_evolution_of_of_op_ot.jpg", width = 10, height = 5) 

ven_mun <- readRDS("ven_mun_elec_2006_2024_final.rds")

view(ven_mun %>% filter(year == 2013) %>% arrange(turnout) %>% head())


vote_rural_urban <- ven_mun %>%
  mutate(op_p = ifelse(year == 2018, NA, op_p)) %>%  # Remove opposition votes in 2018
  group_by(urb_level_ocde_cat) %>% 
  summarise(
    n = n(),
    vote_of = weighted.mean(of_p, validos, na.rm = TRUE),
    vote_op = weighted.mean(op_p, validos, na.rm = TRUE),
    se_vote_of = sd(of_p, na.rm = TRUE) / sqrt(n()),
    se_vote_op = sd(op_p, na.rm = TRUE) / sqrt(n())
  ) %>% 
  drop_na(urb_level_ocde_cat) %>% 
  mutate(urb_level_ocde_cat = factor(urb_level_ocde_cat, levels = c("rural", "intermediate", "urban"), 
                                labels = c("Rural", "Intermediate", "Urban"))) %>%
  pivot_longer(cols = c(vote_of, vote_op), names_to = "Political Block", values_to = "Vote Share") %>%
  mutate(`Political Block` = recode(`Political Block`, "vote_of" = "Chavismo", "vote_op" = "Opposition")) %>%
  ggplot(aes(x = urb_level_ocde_cat, y = `Vote Share`, color = `Political Block`, group = `Political Block`)) + 
  geom_point(position = position_dodge(width = 0.3), size = 3) +  # Slight separation
  geom_errorbar(aes(ymin = `Vote Share` - 1.96 * ifelse(`Political Block` == "Chavismo", se_vote_of, se_vote_op),
                    ymax = `Vote Share` + 1.96 * ifelse(`Political Block` == "Chavismo", se_vote_of, se_vote_op)), 
                width = 0.2, position = position_dodge(width = 0.3)) +  # Match point separation
  scale_color_manual(values = c("Chavismo" = "red", "Opposition" = "blue")) +
  labs(
    title = "Figure 4: Average share of vote by urban-rural classification",
    subtitle = "Venezuela Presidential elections 2006-2024",
    x = "",
    y = "Vote Share (%)",
    color = "Political Bloc",
    caption = "Note: Bars represent 95% confidence intervals.\nThe opposition vote share excludes 2018 due to non-participation.\nClassification based on OECD urban classification: Cities (>50,000 inhabitants, density >1,500/km²), \nIntermediate (>5,000 inhabitants, density >300/km²), Rural areas (low-density or uninhabited)."
  ) + 
  theme_bw() +
  theme(plot.caption = element_text(hjust = 0, face = "italic"),
        axis.title.y = element_text(vjust = 2, size = 12),
                axis.text.x = element_text(size = 10),
                legend.text = element_text(size = 10),
                title = element_text(size = 14))


vote_rural_urban

ggsave(vote_rural_urban, file = "voteshare_ruralurban.jpg", width = 10, height = 6) 


ggarrange(vote_share, vote_rural_urban, nrow = 2, heights = c(2,2)) %>% 
  annotate_figure(., top = text_grob("Figure 2. Vote Share Trends and Electoral Behavior by Urban-Rural Classification in \nVenezuela 2006-2024",
                                     face =  "bold", size = 14)) %>% 
  ggexport(filename = "plots_voteshare_ruralurban.jpg", width = 650, height = 750)

## Map plots ####

ven_elec_mapa <- readRDS("ven_elec_mun_cart.rds")

## Create a separate layer for Esequibo Territory (to overlay hatching lines) ####
esequibo_layer <- ven_elec_mapa %>% filter(entidad == "TERRITORIO ESEQUIBO") %>% 
  mutate(esequibo = "Reclamation")

names(ven_elec_mapa)


## Strongholds

library(RColorBrewer)

ven_elec_mapa <- ven_elec_mapa %>%
  rowwise() %>%
  mutate(
    # Compute mean and median vote share for each party across elections
    mean_of_p = mean(c(of_p_2006, of_p_2012, of_p_2013, of_p_2018, of_p_2024), na.rm = TRUE),
    mean_op_p = mean(c(op_p_2006, op_p_2012, op_p_2013, op_p_2024), na.rm = TRUE),
    mean_otro_p = mean(c(otro_p_2006, otro_p_2012, otro_p_2013, otro_p_2018, otro_p_2024), na.rm = TRUE),
    mean_turnout = mean(c(turnout_2006, turnout_2012, turnout_2013, turnout_2018, turnout_2024), na.rm = T),
    
    median_of_p = median(c(of_p_2006, of_p_2012, of_p_2013, of_p_2018, of_p_2024), na.rm = TRUE),
    median_op_p = median(c(op_p_2006, op_p_2012, op_p_2013, op_p_2024), na.rm = TRUE),
    median_otro_p = median(c(otro_p_2006, otro_p_2012, otro_p_2013, otro_p_2018, otro_p_2024), na.rm = TRUE),
    
    # Stronghold classification based on historical dominance (mean vote share)
    stronghold = case_when(
      mean_of_p > 55 ~ "Chavismo Stronghold",
      mean_op_p > 55 ~ "Opposition Stronghold",
      mean_otro_p > 55 ~ "Other Stronghold",
      TRUE ~ "Battleground"
    )
  ) %>%
  ungroup()

strongholds <- ggplot(ven_elec_mapa) +
  geom_sf(aes(fill = stronghold), color = "black", lwd = 0.1) +
  # Overlay dashed lines specifically over Esequibo territory
  geom_sf(data = esequibo_layer, aes(fill = esequibo), linetype = "dashed", size = 1) +
  scale_fill_manual(
    name = "Political Stronghold",
    breaks = c("Chavismo Stronghold",
               "Opposition Stronghold",
               "Other Stronghold",
               "Battleground"),
    values = c(
      "Chavismo Stronghold" = "red",
      "Opposition Stronghold" = "blue",
      "Other Stronghold" = "purple",
      "Battleground" = "gold"
    )
  ) +
  theme(line = element_blank(),
        axis.text = element_blank(),
        axis.title = element_blank(),
        panel.background = element_blank(),
        plot.caption = element_text(hjust = 0, face = "italic"),
        legend.text = element_text(size = 10),
        title = element_text(size = 14)) +
  labs(
    title = "Figure 2: Strongholds by Municipality",
    subtitle = "Average vote share: Venezuelan presidential elections 2006-2024",
    caption = "Note: Strongholds are municipalities where a party has historically received >55% of the vote.\nMain opposition vote share in 2018 is excluded from the average calculation due to non-participation.\nWhite municipalities correspond to missing data. Grey correspond to the Esequibo territory \nwhich is under reclamation."
  )

strongholds

ggsave("hist_strongholds.jpg", strongholds, width = 10, height = 8)

#### Avergae turnout #####

ven_elec_mapa <- ven_elec_mapa %>% 
  mutate(mean_turnout = case_when(
    is.na(mean_turnout) ~ lag(mean_turnout),
    TRUE ~ mean_turnout
  ))

map_turnout <- ggplot(ven_elec_mapa) + 
  geom_sf(aes(fill = mean_turnout)) +
  geom_sf(data = ven_elec_mapa, fill = NA, color = "black", lwd = 0.15) +
  geom_sf(data = esequibo_layer, fill = "darkgrey", linetype = "dashed", size = 1) +
  theme(line = element_blank(),
        axis.text = element_blank(),
        axis.title = element_blank(),
        panel.background = element_blank(),
        legend.text = element_text(size = 10),
        title = element_text(size = 14)) +
  guides(fill = guide_colorbar(title = "Participation rate (%)")) +  # Colorbar guide for continuous scale
  scale_fill_gradient(
    name = "Turnout \n",
    na.value = "white",
    low = "yellow", 
    high = "darkgreen",
    limits = c(40, 80), # Ensures correct range
    breaks = seq(40, 100, by = 10),  # Adjusts legend steps from 55% to 100%
    labels = function(x) paste0(x, "%")  # Format labels as percentages
  ) +
  labs(
    title = "Figure 1: Average electoral turnout by Municipality",
    subtitle = "Venezuela Presidential elections 2006-2024"
  )

map_turnout

summary(ven_elec_mapa$mean_turnout)

ggsave("map_turnout_2006_2024.jpg", map_turnout, width = 10, height = 8)
